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Cell-centered and node-centered approaches have been compared for unstructured finite- volume discretiza- 
tion of inviscid fluxes. The grids range from regular grids to irregular grids, including mixed-element grids and 
grids with random perturbations of nodes. Accuracy, complexity, and convergence rates of defect-correction 
iterations are studied for eight nominally second-order accurate schemes: two node-centered schemes with 
weighted and unweighted least-squares (LSQ) methods for gradient reconstruction and six cell-centered schemes 
- two node-averaging with and without clipping and four schemes that employ different stencils for LSQ gra- 
dient reconstruction. The cell-centered nearest-neighbor (CC-NN) scheme has the lowest complexity; a version 
of the scheme that involves smart augmentation of the LSQ stencil (CC-SA) has only marginal complexity in- 
crease. All other schemes have larger complexity; complexity of node-centered (NC) schemes are somewhat 
lower than complexity of cell-centered node-averaging (CC-NA) and full-augmentation (CC-FA) schemes. 

On highly anisotropic grids typical of those encountered in grid adaptation, discretization errors of five of 
the six cell-centered schemes converge with second order on all tested grids; the CC-NA scheme with clipping 
degrades solution accuracy to first order. The NC schemes converge with second order on regular and/or 
triangular grids and with first order on perturbed quadrilaterals and mixed-element grids. All schemes may 
produce large relative errors in gradient reconstruction on grids with perturbed nodes. Defect-correction 
iterations for schemes employing weighted least-square gradient reconstruction diverge on perturbed stretched 
grids. Overall, the CC-NN and CC-SA schemes offer the best options of the lowest complexity and second- 
order discretization errors. 

On anisotropic grids over a curved body typical of turbulent flow simulations, the discretization errors 
converge with second order and are small for the CC-NN, CC-SA, and CC-FA schemes on all grids and for 
NC schemes on triangular grids; the discretization errors of the CC-NA scheme without clipping do not con- 
verge on irregular grids. Accurate gradient reconstruction can be achieved by introducing a local approximate 
mapping; without approximate mapping, only the NC scheme with weighted LSQ method provides accurate 
gradients. Defect correction iterations for the CC-NA scheme without clipping diverge; for the NC scheme 
with weighted LSQ method, the iterations either diverge or converge very slowly. The best option in curved 
geometries is the CC-SA scheme that offers low complexity, second-order discretization errors, and fast con- 
vergence. 


I. Introduction 

Both node-centered and cell-centered finite-volume discretization schemes are widely used for complex three- 
dimensional turbulent simulations in aerospace applications. The relative advantages of the two approaches have been 
extensively studied in the search for methods that are accurate, efficient, and robust over the broadest possible range 
of grid and solution parameters. The topic was discussed in a panel session at the 2007 AIAA Computational Fluid 
Dynamics conference, but a consensus did not emerge. One of the difficulties in assessing the two approaches is that 
comparative calculations were not completed in a controlled environment, i.e., computations were made with different 
codes and different degrees of freedom and the exact solutions were not known. 

In this paper, we provide a controlled environment for comparing a subset of the discretization elements needed 
in turbulent simulations, namely that of the inviscid discretization. In particular, we consider a constant-coefficient 
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convection equation as a model for inviscid fluxes. This paper is second in a series of papers on comparison of cell- 
centered and node-centered finite-volume discretizations. It follows Ref. [1], which considered viscous fluxes. The 
ultimate objective of the effort is to construct a uniformly second-order accurate and efficient unstructured-grid solver 
for the Reynolds-Averaged Navier-Stokes equations. 

In this work, we use the method of manufactured solution so that the exact solution is known and conduct com- 
putational studies of accuracy, complexity, and efficiency on two-dimensional grids ranging from structured (regular) 
grids to irregular grids composed of arbitrary mixtures of triangles and quadrilaterals. Highly irregular grids are de- 
liberately constructed through random perturbations of structured grids to bring out the worst possible behavior of the 
solution. Two classes of tests are considered. The first class of tests involves smooth manufactured solutions on both 
isotropic and highly anisotropic grids with discontinuous metrics, typical of those encountered in grid adaptation. The 
second class of tests concerns solutions and grids varying strongly anisotropically over a curved body, typical of those 
encountered in high-Reynolds number turbulent flow simulations. 

There are eight main schemes considered — two representative node-centered schemes with weighted and un- 
weighted least-square methods for gradient reconstruction and six cell-centered schemes. The cell-centered schemes 
include node-averaging schemes with and without clipping and four least-square gradient reconstruction schemes that 
are named according to the stencil used for the least-square fit: a nearest-neighbor scheme uses only face-neighboring 
cells; a smart-augmentation scheme minimally augments the nearest-neighbor stencil; two full augmentation schemes 
with and without weighting use larger stencils that include all node-sharing cells. Each of the schemes considered is 
nominally second-order accurate. 

For the second class of tests, the approximately mapped least-square approach introduced in Ref. [1] is used to 
improve gradient reconstruction accuracy on curved high-aspect-ratio grids. The mapping uses the distance function 
commonly available in practical codes and can be used with any scheme. 

The properties to be compared in this study are computational complexity (operation count) and discretization 
accuracy at equivalent numbers of degrees of freedom as well as convergence rates of defect-correction iterations with 
a first-order driver. The effect of clipping is studied for the node-averaging schemes. 

The material in this paper is presented in the following order. Section II introduces the computational grids used 
in the current study. A brief explanation of finite-volume discretizations in Section III is followed by the estimates of 
discretization complexity for two- and three-dimensional grids given in Section IV. Section V outlines the analysis 
methods used in this study. A brief introduction of the model equation in Section VI precedes results provided in 
Section VII on accuracy of finite-volume solutions and gradients and on convergence rates of defect-correction iter- 
ations observed on isotropic irregular grids. The effect of clipping on accuracy of node-averaging schemes is also 
studied in this section. Section VIII compares the finite-volume discretizations on stretched highly anisotropic grids 
in rectangular geometries. Section IX provides comparisons for irregular high-aspect-ratio grids in curved geometries. 
Conclusions and recommendations are offered in Section X. 

II. Grids 

This paper studies finite-volume discretization (FVD) schemes for inviscid fluxes on grids that are loosely defined 
as irregular. A grid is classified as regular if it can be derived by a smooth mapping from a grid with (1) a periodic 
node connectivity pattern (i.e., the number of edges per node changes periodically) and (2) a periodic cell distribution 
(i.e., the grid is composed of periodically repeated combinations of cells). Regular grids include, but are not limited to, 
grids derived from Cartesian ones - triangular grids obtained by diagonal splitting with a periodic pattern, smoothly 
stretched grids, skewed grids, smooth curvilinear grids, etc. Grids that are not regular are called irregular grids. We 
are especially interested in unstructured grids, e.g., grids with the number of edges changing from node to node with 
no pattern. 

The regular and irregular grids considered in this paper are derived from an underlying (possibly mapped) Cartesian 
grid with mesh sizes h x and h y and the aspect ratio A = h x /h y ; both mesh sizes of the underlying grid are assumed to 
be small, h y 1 , h x <C 1 . Irregularities are introduced locally and do not affect grid topology and metrics outside of a 
few neighboring cells. A local grid perturbation is called random if it is independent of local perturbations introduced 
beyond some immediate neighborhood. For computational grids generated for the reported studies, local and random 
grid irregularities are introduced in two ways: (1) the quadrilateral cells of the underlying grid are randomly split (or 
not split) into triangles; (2) the grid nodes are perturbed from their original positions by random shifts, where the shifts 
are fractions of a local mesh size. 

Four basic grid types are considered: (I) regular quadrilateral (i.e., mapped Cartesian) grids; {II) regular tri- 
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angular grids derived from the regular quadrilateral grids by the same diagonal splitting of each quadrilateral; (III) 
random triangular grids , in which regular quadrilateral are split by randomly chosen diagonals, each diagonal orienta- 
tion occurring with probability of half; (IV) random mixed-element grids , in which regular quadrilateral are randomly 
split or not split by diagonals; the splitting probability is half; in case of splitting, each diagonal orientation is chosen 
with probability of half. Nodes of any basic-type grid can be perturbed from their initial positions by random shifts, 
thus leading to four additional perturbed grid types which are designated by subscript p as (I p )-(IV p ). Grids of types 
(III) — (IV) and (II Ip) — (IV P ) are irregular (and unstructured) because there is no periodic connectivity pattern. 
All perturbed grids are irregular because there is no periodic cell distribution. The representative grids are shown in 
Figure 1. 


(a) Type (I): regular quadri- 
lateral grid. 



(e) Type (Ip): perturbed 

quadrilateral grid. 



(b) Type (II): regular trian- 
gular grid. 



(f) Type (II P ): perturbed 
triangular grid. 



(c) Type (III): random tri- 
angular grid. 



(g) Type (I I I p ): perturbed 
random triangular grid. 



(d) Type (IV): random 

mixed grid. 



(h) Type (IV P ): perturbed 
random mixed grid. 


Figure 1. Typical regular and irregular grids. 

Our main interest is the accuracy of FVD schemes on general irregular (mostly unstructured) grids with a minimum 
set of constraints. In particular, we do not require any grid smoothness, neither on individual grids nor in the limit 
of grid refinement. The only major requirement for a sequence of refined grids is to satisfy the consistent refinement 
property. The property requires the maximum distance across the grid cells to decrease consistently with increase 
of the total number of grid points, N. In particular, the maximum distance should tend to zero as iV -1 / 2 in 2D 
computations. For 3D unstructured grids, the consistent refinement property has been studied elsewhere. 2 On 2D 
grids, the effective mesh size , h e , is computed as the L\ norm of the square root of the control volumes. 

The locations of discrete solutions are called data points. For consistency with the 3D terminology, the 2D cell 
boundaries are called faces, and the term “edge” refers to a line, possibly virtual, connecting the neighboring data 
points. Each face is characterized by the directed-area vector , which is directed outwardly normal to the face with the 
amplitude equal to the face area. 

The random node perturbation in each dimension is defined as jph , where p € [—1, 1] is a random number and h 
is the local mesh size along the given dimension. With these perturbations, triangular cells in the rectangular geometry 
can approach zero volume. The random perturbations are introduced independently on all grids in grid refinement 
implying that grids of types (I p ) — ( IV P ) are grids with discontinuous metrics, e.g., ratios of neighboring cell volumes 
and face areas are random on all grids and do not approach unity in the limit of grid refinement. 


III. Finite- volume discretization schemes 

The FVD schemes are derived from the integral form of a conservation law 
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( 1 ) 


on n 

where Q is a control volume, F is the flux through the boundary dfl, n is the outward unit normal vector, and / is 
a force function. The general FVD approach requires partitioning the domain into a set of non-overlapping control 
volumes and numerically implementing equation (1) over each control volume. 



Figure 2. Control-volume partitioning for finite- volume discretizations. Numbers 0 — 12 and letters A — L denote grid nodes and primal 
cell centers, respectively. The control volume for a node-centered discretization around the grid node 0 is shaded. The control volume for 
a cell-centered discretization around the cell center A is hashed. 

Cell-centered (CC) discretizations assume solutions are defined at the centers of the primal grid cells with the 
primal cells serving as the control volumes. The cell center coordinates are typically defined as the averages of 
the coordinates of the cell’s vertexes. Note that for mixed-element grids cell centers are not necessarily centroids. 
Node-centered (NC) discretizations assume solutions are defined at the primal mesh nodes. For NC schemes, control 
volumes are constructed around the mesh nodes by the median-dual partition: the centers of primal cells are connected 
with the midpoints of the surrounding faces. These non-overlapping control volumes cover the entire computational 
domain and compose a mesh that is dual to the primal mesh. Both cell-centered and node-centered control-volume 
partitions are illustrated in Figure 2. 

The fluxes at a control-volume face are computed according to the Roe scheme, 3 

(F • n) = i [(F r ■ n) + (F l • n)] - ± |A| (Qr - Ql), (2) 

where, Ql and Qr are the “left” and “right” solution reconstructions; Fl and Fr are the corresponding “left” and 
“right” numerical fluxes; | A| is the Roe’s approximate Riemann solver matrix. The solutions Ql and Qr are linearly 
reconstructed at the face by using solutions defined at the control volume centers and solution gradients reconstructed 
at each control volume. Various FVD schemes differ in the way they reconstruct gradients at the control volumes. 

For cell-centered schemes, the face-based flux integration over a control-volume face is approximated as the inner 
product of F computed at the face center and the face directed area vector. The integration scheme is second-order 
accurate on grids of all types. For node-centered schemes, the edge-based flux integration scheme approximates the 
integrated flux through the two faces linked at an edge midpoint by multiplying F computed at the edge midpoint 
with the combined-directed-area vector, n = riL + iir, where riL and iir are directed-area vectors of the left and 
right faces, respectively. The integration scheme is computationally efficient and second-order accurate on regular 
and triangular grids of types (I), (II), (HI), (Up)* and (III p ); the integration accuracy degenerates to first order on 
mixed-element and perturbed quadrilateral grids of types (IV), (IV P ), and (I p ): 2,4,5 

The forcing term integration over the control volume is approximated as the value at the control-volume center 
multiplied by the volume Q|. This approximation is second-order accurate when the control-volume center coincides 
with the centroid. On general irregular grids, the control-volume center is not necessarily the centroid, and the approx- 
imation becomes locally first-order accurate. However, with grid irregularities introduced locally and randomly (thus. 
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implying a zero-mean distribution of the deviations between control-volume centers and centroids), the integral of the 
forcing term over any sub-domain of size 0(1) is approximated with second order. 


A. Cell-centered schemes 

1. Node averaging schemes 

In the cell-centered node-averaging (CC-NA) schemes, the solution values are first reconstructed at the nodes from the 
surrounding cell centers. With respect to Figure 2, the solution at the node 0 is reconstructed by averaging solutions 
defined at the cell centers A, B , and C. The solution reconstruction proposed in Refs. [6, 7] and used in Ref. [8] is an 
averaging procedure that is based on a constrained optimization to satisfy some Laplacian properties. The scheme is 
second-order accurate and stable when the coefficients of the introduced pseudo-Laplacian operator are close to 1 . It 
has been shown 9 that this averaging procedure is equivalent to an unweighted least-square linear fit. 

The gradient at the cell Q is reconstructed by the Green-Gauss formula, 



on 


Uhds, 


( 3 ) 


where |f2| is the cell volume, n is the outward unit normal, ds is the area differential, and integration is performed 
over the cell boundary, Oil. For second-order accuracy, the solution at a face is computed by averaging the values at 
the face nodes and the integral over the face is approximated by the product of the solution and the face directed area. 

On highly stretched and deformed grids, some coefficients of the pseudo-Laplacian may become negative or larger 
than 2, which has a detrimental effect on stability and robustness. 10- 11 Holmes and Connell 6 proposed to enforce 
stability by clipping the coefficients between 0 and 2. The CC-NA schemes with clipping (CC-NA-CLIP) represent 
a current standard in practical computational fluid dynamics for applications involving cell-centered finite volume 
formulations. 12 As shown further in the paper, clipping seriously degrades accuracy of the solutions and gradients. 


2. Least-square schemes 

An alternative approach relies on a least-square method for gradient reconstruction, in which the linear approximation 
obtained at a control volume is required to coincide with the solution value at the control volume center. In this paper, 
both weighted and unweighted least-square methods are considered. The weighted method is designated as WLSQ 
herein and the unweighted method is used as default without designation. In the WLSQ method, the contributions to 
the minimized functional are weighted with weights inversely proportional to the distance from the control-volume 
center. In the unweighted method, all contributions are equally weighted. 

The stencils used in the gradient fits are discussed with respect to Figure 2. Three types of stencils are considered 
— nearest neighbor (NN), full augmentation (FA), and smart augmentation (SA) stencils. The NN stencil involves 
only centers of face-neighbor cells; the FA stencil includes all the cells that share a vertex with the given cell, i.e., all 
the cells involved in CC-NA gradient reconstruction; the SA stencil is an adaptive stencil that provides a minimally 
necessary extension of the NN stencil to improve convergence rates of the defect-correction iterations (DCI) with the 
first-order cell-centered FVD scheme as the driver. For cell-center A , the NN stencil includes neighbors B 1 C, D, 
and E\ the FA stencil includes additionally neighbors F, G , II . /, J, K, and L; the SA applies an augmentation test to 
the NN stencil and expands it only if necessary and by choosing only appropriate cells from the augmentation pool 
provided by the FA method. 

Initially, the CC-SA scheme is identical to the CC-NN scheme. In stencil augmentation at each cell, the augmen- 
tation test computes the quantity Ci c = 1 1 — dsA/di |, where dsA and d\ are the respective main-diagonal coefficients 
of full linearizations of the current CC-SA and the first-order driver schemes for a constant-coefficient convection 
operator. The test is applied for a preselected number of representative convection directions indexed by i c . In the 
algorithm implemented for this paper, the current CC-SA scheme is considered sufficiently augmented if the augmen- 
tation indicator 


AI = max C lc < e, (4) 

where e = 0.4 is a user-defined tolerance. Smaller values of e imply larger CC-SA stencils. If augmentation is 
required, only one cell from the augmentation pool is added to the stencil. The cells from the pool are tested one by 
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one until a cell that brings AI below the e-threshold is found. If no such single cell has been found, the cell that makes 
the best improvements in AI is added to the stencil, and the augmentation procedure repeats. Note that it is possible 
that at the end, the user-defined tolerance has not been achieved. Even in these instances, the smart augmentation adds 
only cells that reduce AI, thus, providing a much smaller stencil than CC-FA stencil even in the worst-case scenario. 
Note, also, that the results of smart augmentation may depend on the order in which cells have been augmented. In the 
current paper, a sequential smart augmentation order has been used, while a fully parallel version which is independent 
on the augmentation order has also been developed and implemented. 

B. Node-centered schemes 

For the node-centered computations, the current standard employs a least-square gradient reconstruction. The typical 
stencil at a control volume involves all nodes linked by an edge. For example, with reference to Figure 2, the least- 
square fit for the shaded control volume centered at node 0 includes nodes 1,2, and 4. Both weighted and unweighted 
least-square methods are evaluated. 


IV. Complexity 


A. Flux integration complexity 

In this section, the complexity associated with flux integration in 3D cell-centered or node-centered FVD schemes 
is estimated. The complexity is measured as the number of flux-reconstruction instances required for one residual 
evaluation. Flux reconstructions are the main contributers to the operation counts associated with flux integration; other 
aspects of the discretization, such as determining the solution values or solution-gradient values require additional 
considerations. Three types of primal meshes are considered: (1) fully-tetrahedral, (2) fully-prismatic, (3) fully- 
hexahedral. 

An underlying Cartesian grid is considered and split into the various elements. The splitting into tetrahedra assumes 
each hexahedral defined by the grid is split into 5 tetrahedra with one of the tetrahedra being completely interior to the 
hexahedral (i.e., its faces are not aligned with any of the hexahedral faces - see Figure 3). Note that there are other 
partition strategies that lead to different number of tetrahedra per hexahedral; for example, dividing the hexahedral 
into two triangular prisms with subsequent division of each of the prisms into 3 tetrahedra leads to 6 tetrahedra per 
hexahedral. In this section we do not consider other possible partitions. 



Figure 3. Splitting hexahedral into 5 tetrahedra. 

Table 1 shows complexity estimates for two node-centered and one cell-centered 3D FVD schemes. Only interior 
discretizations are estimated; boundary effects are neglected. Both node-centered discretizations assume a median- 
dual partition of the domain. In such a partition, the constituent dual control volumes are bounded by generally 
non-planar dual faces formed by connecting 3 types of points: (1) edge midpoints, (2) element-face centroids, and 
(3) element centroids. FVD schemes with edge-based flux integration, such as NC schemes used in the current study, 
approximate integration over all of the constituent dual faces surrounding an edge midpoint by evaluating the flux 
at the edge midpoint; the directed area is taken as the combined directed area. FVD schemes with face-based flux 
integration reconstruct fluxes at each of the constituent dual faces separately and use local directed areas. For the 
present estimation, we assume that each flux-reconstruction instance requires the same operation count, in particular, 
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the approximate Riemann solver is applied at each reconstruction point. In fact, significant savings can be achieved, if 
the dissipation matrix is computed once for all control surfaces surrounding an edge. The first node-centered scheme 
is a linear 3D FVD scheme with edge-based flux integration; the second node-centered scheme is a linear 3D FVD 
scheme with face-based flux integration. The cell-centered formulation uses a face-based flux integration scheme with 
one flux reconstruction per control face. 

Two estimates of complexity are given. The first estimate assumes that any constituent quadrilateral face in the 
control surface is broken into two triangular faces. The second estimate (in parentheses) assumes any constituent 
quadrilateral face is approximated as planar. The former is required to ensure a precise (water-tight) definition of the 
control surface and can serve as a measure of the complexity in integration of the physical flux terms. The latter can 
serve as an estimate of the complexity associated with numerical dissipation terms, in which details of the control- 
surface can be neglected. 


Elements 

Cell -centered 
face-based flux integration 

Node-centered 
edge-based flux integration 

Node-centered 
face-based flux integration 

Tetrahedral 

4(4) 

12 

120 (60) 

Prismatic 

8(5) 

8 

72 (36) 

Hexahedral 

12 (6) 

6 

48 (24) 


Table 1. Number of ftix-reconstruction instances per equation for 3D FVD discretizations. 

The complexities of cell-centered and node-centered FVD schemes with edge-based flux integration are reasonably 
close. Unfortunately, as shown in this paper and also previously, 2- 4-5 the accuracy of the edge-reconstruction FVD 
scheme degenerates to first order on perturbed quadrilateral and general mixed-element grids. To maintain the second- 
order accuracy on general grids, one can employ the node-centered scheme with face-based flux integration, but the 
integration complexity of this formulation substantially exceeds the complexity of the cell-centered FVD scheme. 
These results are in agreement with the observations made by Delanaye and Liu 13 leading to the selection of a cell- 
centered discretization. 

B. Size of inviscid stencil 

Another important measure of complexity of an FVD scheme is the size of the full-linearization stencil. The size of 
the 2D and 3D full-linearization stencil is examined for the inviscid cell-centered and node-centered FVD schemes. 
Cartesian meshes are split into triangular and tetrahedral elements, as in the previous section, again neglecting bound- 
ary effects. Estimates are compared to numerical calculations on an actual 3-D grid that includes boundary effects; the 
grid is a viscous fully-tetrahedral grid composed of 16,391 nodes. 

In three dimensions, half of the grid nodes have 18 adjacent edges (32 adjacent tetrahedra) and half have 6 adjacent 
edges (8 adjacent tetrahedra). Each of the tetrahedra interior to an originally-hexahedral cell is defined by four nodes, 
each with 18 adjacent edges. Each of the four surrounding tetrahedra within an originally-hexahedral cell is defined 
by three nodes with 18 adjacent edges and 1 node with 6 adjacent edges. 

For reference. Table 2 shows the average and maximum number of edges, n e d ge , connecting to a grid node. The 
average number of connecting edges sets the least-square stencil size for the node-centered scheme as n e d ge + 1. The 
number of connecting edges is also an important factor for the CC-NA schemes because it characterizes the number 
of elements sharing the node and therefore the number of cells used for averaging data to the grid node. Generally 
speaking, the number of edges is not bounded in 3D and, thus, the corresponding CC-NA stencil size is not bounded. 


Dimension 

nedge (Average) 

nedge (Maximum) 

2D 

6 

8 

3D 

12 

18 


Table 2. Edges connecting to a grid node in the split Cartesian grids. 

For the inviscid discretization, the DCI with a first-order driver is generally used to converge the residual; thus, 
it is important to consider first-order and second-order linearizations. For the first-order cell-centered FVD scheme. 
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the size of the linearization stencil is simply the number of faces plus one (to account for the central node). For the 
first-order node-centered discretization, the size of the linearization stencil is the number of edges connecting to a node 
plus one. Table 3 shows 2D and 3D linearization stencil sizes. The cell-centered discretization has nearly a factor of 3 
smaller stencil in 3D. 


Elements 

Node-centered 

Cell-centered 

Estimate 2D 

7 

4 

Estimate 3D 

13 

5 

Numerical 3D 

14 

5 


Table 3. Average size of the inviscid first-order FVD stencil on triangular/tetrahedral grids in 2D/3D. 

For second-order accuracy, all schemes reconstruct gradients at the control volumes. The node-centered discretiza- 
tions use a least-squares approach and require solutions at the neighbor-of-neighbor nodes and a correspondingly large 
linearization stencil. The cell-centered CC-NA schemes have even larger linearization stencils which include all cells 
contributing to solution reconstruction at any node of a face-neighboring cell. Stencils of CC-FA schemes are the same 
as CC-NA stencils. The CC-NN stencil also uses a least-squares approach to fitting the gradient in reconstruction, but 
requires a much smaller stencil which includes only neighbor-of-neighbor cells. Table 4 shows stencil sizes for 2D 
and 3D; in 3D, only the splitting shown in Figure 4 is considered. In three-dimensions, the NC stencil is significantly 
smaller than the CC-NA and CC-FA stencils. In both 2D and 3D, the CC-NN stencil is the smallest. 


Elements 

NC 

CC-NA 

CC-NN 

Estimate 2D 

23 

25 

9 

Estimate 3D 

75 

139 

15 

Numerical 3D 

63 

118 

15 


Table 4. Average size of the inviscid second-order stencil for 2D/3D discretizations with triangular/tetrahedral elements. 


The numbers are so striking that it is useful to show the stencils for a single shaded control volume in Figure 4 for 
each approach. The stencil sizes are 25, 25, and 9 for the NC, CC-NA, and CC-NN schemes, respectively. Note that 
the stencil size for the NC control-volume adjacent to the one shown in Figure 4 is 21; thus, the average of 23 is shown 
in Table 4. Also, for the 3D NC schemes, the nodes with 6 and 18 edges have stencil sizes of 57 and 93, respectively; 
thus, the average of 75 is shown in the table. For the CC-NA and CC-FA schemes, the cells at the corners of the 
original Cartesian cell have a stencil size of 149 and those fully interior to the original Cartesian cell have a stencil 
size of 99. Since there is one interior tetrahedron for each of the four corner tetrahedrons, the average of 139 is shown 
in the table. 



(a) NC scheme. 




(b) CC-FA and CC-NA schemes. 


(c) CC-NN scheme. 


Figure 4. Inviscid 2D stencil for shaded control volume. 


8 of 24 


American Institute of Aeronautics and Astronautics 








V. Analysis 


A. Method of manufactured solution 

Accuracy of FVD schemes is analyzed for known exact or manufactured solutions. The forcing function and boundary 
values are found by substituting this solution into the governing equations, including boundary conditions. The discrete 
forcing function is defined at the data points. 


1. Discretization error 


The main accuracy measure is the discretization error, E d , which is defined as the difference between the exact discrete 
solution, U h , of the discretized equations (1) and the exact continuous solution, U, to the corresponding differential 
equations 

E d = U~ U h ; (5) 


U is sampled at data points. 


2. Truncation error 

Another accuracy measure commonly used in computations is truncation error. Truncation error, E t , characterizes the 
local accuracy of approximating the differential equations. For finite differences, it is defined as the residual obtained 
after substituting the exact solution U into the discretized differential equations. 14 For FVD schemes, the traditional 
truncation error is usually defined from the time-dependent standpoint. 15 - 16 In the steady-state limit, it is defined (e.g., 
in Ref. [17]) as the residual computed after substituting U into the normalized discrete equations (1), 



- Jf f h dn + <j> (F h ■ n) ds , 

n on 


where |0| is the measure of the control volume. 


(6) 


|fi| = JJ dfl, (7) 

n 

F h is a numerical flux evaluated at the control-volume boundary (El, f h is an approximation of the forcing function 
f on O, and the integrals are computed according to some quadrature formulas. Note that convergence of truncation 
errors is expected to show the order property only on regular grids; on irregular grids, it has been long known that 
the design-order discretization-error convergence can be achieved even when truncation errors exhibit a lower-order 
convergence or, in some cases, do not converge at all. 17-21 


3. Accuracy of gradient reconstruction 


Yet another important accuracy measure is the accuracy of gradient approximation at a control-volume. For second- 
order convergence of discretization errors, the gradient is usually required to be approximated with at least first order. 
For each control- volume, accuracy of the gradient is evaluated by comparing the reconstructed gradient, V r , with the 
exact gradient, V exa ct> computed at the control-volume center. The accuracy of gradient reconstruction is measured 
as the relative gradient error: 


E re 1 — 


lie'll’ 


( 8 ) 


where functions e and G are amplitudes of the gradient error and the exact gradient, respectively, evaluated at face 
centers; 


e = | V r U h - V exact U 1 1 and G = |V exac t£7|; 


(9) 


U and U h are a differentiable manufactured solution and its discrete representation (usually injection) on a given grid, 
respectively; || • || is a norm of interest computed over the entire computational domain. 
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4. Convergence of iterative solvers 


Besides accuracy, an important quality of a practical discretization is availability of an affordable solver. For FVD 
schemes with low complexity, such as CC-NN and CC-SA, an efficient solution method would use a full linearization 
in relaxation of the target FVD scheme. For FVD schemes with high complexity, such as CC-NA, CC-FA, and even 
NC schemes, iterations with the full linearization are not affordable; DCI schemes with linearized first-order drivers 
are common methods used in practical computations. In this view, stability and convergence rates of DCI are also 
analyzed. Let u h be the current solution approximation. The DCI method is defined in the following two steps: 

1. The correction v h is calculated from 

L h d v h = R h (u h ) , (10) 

where R h ( u h ) is the residual of the target FVD scheme and L d is a driver scheme. 

2. The current approximation is corrected 

u h = u h + v h . (11) 

All considered second-order FVD schemes use the first-order upwind FVD scheme as a driver. 


VI. Convection equation 

The linear convection equation 

(a • V) (7 = / (12) 

is considered as a model for inviscid fluxes; a is a vector-function of spatial variables. The forcing function / is 
independent of the solution U. Boundary conditions are typically defined either in a weak form as the normal flux, 
(F • n) = U (a • n), given at the inflow boundary or as over-specified conditions, in which solutions at control volumes 
that include nodes edge-connected to the boundary are over-specified from the manufactured solution. In the tests 
reported further in this paper, the convection direction is constant, a = (sin (yg) , cos (jg)), and boundary conditions 
are over-specified. 


VII. Isotropic irregular grids 


A. Grid refinement 

All computations in this section are performed for for the manufactured solution U = — cos (27ra: — iry). Sequences 
of consistently refined grids of types ( III p ) and ( IV v ) are generated on the unit square [0, 1] x [0, 1]. Irregularities 
are introduced at each grid independently, so the grid metrics remain discontinuous on all the grids. The ratio of 
areas of neighboring faces can be as large as 3-\/2; because a control volume can be arbitrarily small, the ratio of the 
neighboring volumes can be arbitrarily high. Two node-centered and six cell-centered schemes are considered; NC, 
NC-WLSQ, CC-SA, CC-NN, CC-FA, CC-FA- WLSQ, CC-NA and CC-NA-CLIP. On grids of type ( III p ), CC-SA 
scheme augments about 50% of the interior least-square stencils and CC-NA-CLIP clips about 10% of the interior 
nodes. On grids of type ( IV P ), CC-SA scheme augments between 25% and 30% of the interior least-square stencils 
and CC-NA-CLIP clips about 3% of the interior nodes. On grids of both types, about 80% of the augmented stencils 
increase the stencil size just by one cell, about 20% by 2 cells, and less than 1% by more than 2 cells. 

B. Gradient reconstruction accuracy 

For second-order discretization accuracy, the gradient reconstruction is required to be at least first-order accurate. 
To evaluate the gradient reconstruction accuracy, the computational gradients have been reconstructed within interior 
control volumes from the manufactured solution evaluated at the data points and compared with the exact gradients 
computed at the control-volume centers. Figure 5 shows convergence of the L a c norms of relative gradient errors on 
grids of types (I Up) and ( IV P ). Only errors computed with the CC-NA-CLIP scheme do not converge in grid refine- 
ment. Similar absence of convergence has been observed and reported previously 1 for gradients reconstructed with the 
clipped CC-NA scheme within control-volume faces. All other methods provide first-order gradient approximations 
on grids of both types. 
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(a) Grids of type ( III V ). 


(b) Grids of type ( IV P ). 


Figure 5. Accuracy of gradient reconstruction for cell-centered FVD schemes on isotropic irregular grids. Manufactured solution is 

U = — cos ( 2irx — try). 


C. Convergence of truncation and discretization error 

Numerical tests evaluating convergence of truncation and discretization errors are performed for the constant-coefficient 
convection equation (12). Figures 6 and 7 show convergence of the L \ norms of truncation and discretization errors, 
respectively. 

Truncation errors of all the cell-centered schemes (except the CC-NA-CLIP scheme ) converge with first order on 
grids of both types and truncation errors of the node-centered schemes converge with first order on triangular grids of 
type (///p); the corresponding discretization errors converge with second order. As predicted in Refs. [2,5], truncation 
errors of node-centered schemes do not converge on mixed-element grids; discretization errors converge with first 
order. The reason for this convergence degradation is the edge-based flux integration scheme, which is second-order 
accurate on simplex (triangular and tetrahedral) grids, but only first-order accurate on perturbed quadrilateral and 
general mixed-element grids. As shown in Ref. [5], with a more accurate face-based flux integration scheme, second- 
order accuracy is achieved with NC schemes on arbitrary grids. Although barely discernible, convergence of truncation 
and discretization errors of the CC-NA-CLIP scheme deteriorates on finer grids. Detailed tests performed on finer 
grids and reported in a subsequent section show that truncation error convergence stagnates and discretization error 
convergence deteriorates to first order. Also not shown, convergence of the L 0 0 norms of the CC-NA-CLIP scheme 
show signs of deterioration on coarser grids. For other schemes, convergence slopes are the same for all norms and do 
not change on finer grids. 

All second-order discretization error plots are very close to each other indicating similar accuracy on grids with 
equivalent number of degrees of freedom. For reference. Figures 7(a) and 7(b) include the convergence plots of 
“ideal” discretization errors computed with the CC-EG scheme that uses exact gradients evaluated at each cell from 
the manufactured solution. These plots represent the best-possible second-order convergence, which can be achieved 
on given grids. Close proximity of the actual and the ideal second-order discretization errors indicates that the accuracy 
is nearly optimal. 

D. Convergence of defect-correction iterations 

Convergence of DCI is studied for the second-order FVD schemes on isotropic grids of types ( III p ) and ( / V p ) with 
65 2 nodes. The forcing term and the boundary conditions are set to zero. The initial solution is random. Convergence 
rates are shown in Figure 8. As was mentioned above, the CC-SA and CC-NN schemes have small stencils and can 
be relaxed with full linearization of target second-order operators. However for consistency, convergence rates of DCI 
are shown for these schemes as well. 

The DCI method for all schemes converges fast with an average convergence rate per iteration better than 0.6. The 
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Figure 6. Convergence of Li -norms of truncation errors of FVD schemes on irregular grids. Manufactured solution is U = 

— cos (2 ■kx — iry). 
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(a) Grids of type (III P ). 


(b) Grids of type ( IV P ) 


Figure 7. Convergence of L \ -norms of discretization errors of FVD schemes on irregular grids of of types (III P ) and (I V p ). Manufactured 
solution is U = — cos (2nx — ny). 


convergence plots can be divided into three parts: initial convergence, transition, and asymptotic convergence. Initial 
convergence is typically fast for random initial solutions. The number of iterations transitions within the transition 
region grows slightly on finer grids. Asymptotic convergence rates for all schemes are around 0.5 per iteration. Note, 
that on grids of type (/), all studied discretization schemes correspond to the Fromm discretization of the convection 
equation. A detailed study of DCI for the Fromm discretization on Cartesian grids has been reported elsewhere. 22 
Note, also, that reported problems with stability of DCI for the WLSQ schemes 23 and for the CC-NA scheme without 
clipping 6 are not evident on these isotropic grids. 

E. Effects of clipping 

The tests reported in this section are performed for the CC-NA and CC-NA-CLIP schemes and demonstrate detrimental 
effects of clipping on convergence of gradient-reconstruction, truncation, and discretization errors in grid refinement. 
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(a) Grids of type ( III P ). 


(b) Grids of type ( IV P ) 


Figure 8. Convergence of Lj -norms of residuals in DCI for second-order FVD schemes with first-order drivers on isotropic irregular grids 
of types (I I I p ) and (IV P ). 


Considered irregular triangular grids of type {IIIp) are characterized by a higher percentage of clipped nodes; about 
10% of the interior nodes are clipped. Figure 9(a) shows an example of a grid of type ( III p ) with 17 2 nodes; nodes 
where clipping occurs are circled. 

Figure 9(b) shows that the gradients reconstructed by the CC-NA-CLIP scheme do not approximate the exact gra- 
dients. The CC-NA scheme provides a first-order accurate gradient reconstruction, which is sufficient for second-order 
discretization accuracy. Figures 9(c) and 9(d) exhibit convergence of the L\ norms of truncation and discretization 
errors, respectively. The CC-NA scheme demonstrates first-order convergence of truncation errors and second-order 
convergence of discretization errors. Truncation errors are very similar on coarse grids, but start to diverge on finer 
grids. Truncations errors of the CC-NA scheme demonstrate clear first-order convergence; truncation errors of the CC- 
NA-CLIP scheme converge slower on finer grids and eventually stagnate. The discretization error convergence of the 
CC-NA-CLIP scheme exibits second order on the coarse grids, but then degrades to first order. Although not shown, 
the Loo norm of discretization errors of the CC-NA-CLIP scheme shows degradation on coarser grids in grid refine- 
ment; asymptotically, L a a norms of both node-averaging schemes converge with the same orders as the corresponding 
L i norms. Note that on grids with a small percentage of clipped nodes, convergence degradation becomes visible only 
on very fine grids. This may explain why such degradation has not been reported for practical computations. 

VIII. Anisotropic irregular grids 


A. Grid stretching 

In this section, we study FVD schemes on stretched grids generated on rectangular domains. Figure 10 shows an 
example grid of type ( III p ) with the maximal aspect ratio A = 10 3 . The manufactured solution is U = sin(7riE-p27r y). 
A sequence of consistently refined stretched grids is generated on the rectangle (x, y) £ [0, 1] x [0, 0.5] in the following 
3 steps. 

1. A background regular rectangular grid with N = (N x + 1) x (N y + 1) nodes and the horizontal mesh spacing 
h x = jj- is stretched toward the horizontal line y = 0.25. The //-coordinates of the horizontal grid lines in the 
top half of the domain are defined as 

y%L + 1=0.25; y 0 =y J - 1 + hyf3 J ~(— +1 ) , j = ^ + 2, . . . , N y , N y + 1. (13) 

Here h y = is the minimal mesh spacing between the vertical lines; A = 10 3 is a fixed maximal aspect ratio; 
/3 is a stretching factor, which is found from the condition t/jv+i = 1- The stretching in the bottom half of the 
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Figure 9. Accuracy of CC-NA schemes on isotropic irregular triangular grids of type ( III P ). Manufactured solution is U = 

— cos (2nx — Try). 


domain is defined analogously. 

2. Irregularities are introduced by random shifts of interior nodes in the vertical and horizontal directions. The 
vertical shift is defined as A yj = ^pmin(/i^ -1 , h J y ), where p is a random number between —1 and 1, and /i^ -1 
and h J y are vertical mesh spacings on the background stretched mesh around the grid node. The horizontal shift 
is introduced analogously, A Xi = \ph x . With these random node perturbations, all perturbed quadrilateral cells 
are convex. 

3. Each perturbed quadrilateral is randomly triangulated with one of the two diagonal choices; each choice occurs 
with a probability of one half. 

B. Gradient reconstruction accuracy 

A recent study 24 assessed accuracy of gradient approximation on various irregular grids with high aspect ratio A = 

jf- 1. The study indicates that for rectangular geometries and functions predominantly varying in the direction of 

rly 
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Figure 10. Random triangular stretched grid with 17 x 65 nodes. 


small mesh spacing i //-direction ), gradient reconstruction is accurate. For manufactured solutions significantly varying 
in the direction of larger mesh spacing (x-direction), the gradient reconstruction may produce extremely large 0(Ah x ) 
relative errors affecting the accuracy of the //-directional gradient component. Figure 1 1 shows examples of first-order 
accurate gradient approximations that exhibit large relative errors on high-aspect-ratio grids of type (III). 
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(a) Aspect ratio A = 10 6 . 


(b) Aspect ratio A = 10 3 . 


Figure 11. Relative errors in approximation of gradients for the manufactured solution U = sin(7i\>: + 2 Try) on anisotropic grids of type 
(III) downscaled toward the focal point ( x , y) = (0.3, 0.5). 

Evaluation of gradient reconstruction accuracy is performed with the methodology of downscaling described in 
detail elsewhere. 2,5 The computational tests are performed on a sequence of downscaled narrow domains L x (L/A) 
centered at the focal point (x,y) = (0.3, 0.5). The scale L changes as L = 2~ n ,n = 0, 8 and the considered 
aspect ratios are A = 10 6 and A = 10 3 ; the latter corresponds to the highest aspect ratio observed at the central 
line of the stretched grid shown in Figure 10. On each domain, an independent high-aspect-ratio random grid of type 
(III) with 9 2 nodes is generated; the grid aspect ratio is fixed as A on all scales. The gradient reconstruction accuracy 
was measured at the interior control volumes. Only weighted-least-square schemes, NC-WLSQ and CC-FA-WLSQ, 
provide accurate gradients, the relative errors of gradient reconstructions provided by all other schemes are several 
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orders of magnitude larger, directly proportional to the aspect ratio A, and converge with first order. 

A summary of the results concerned with gradient accuracy on anisotropic grids is presented in Table 5. All 
considered gradient reconstruction methods may generate large relative errors on perturbed grids of types ( I p ) — (IV P ). 
Only the NC-WLSQ scheme provides gradient reconstruction accuracy on all unperturbed grids. On perturbed grids, 
there are topologies, where all stencil points are almost equidistant from the stencil center, and the WLSQ method is 
ineffective. Such situations occur more frequently for cell-centered schemes; all cell-centered schemes may generate 
large gradient errors even on unperturbed mixed-element grids of type (IV). The CC-NN, CC-NA, and CC-FA- 
unweighted methods may also have large relative errors on random triangular grids of types ( III ); the CC-FA-WLSQ 
method always provides accurate gradients on these grids. 

Table 5. Relative error of gradient reconstruction. 


Grids 

(I) 

(II) 

(III) 

(IV) 

(i P ) - (iv P ) 

NC 

o(hl) 

o(hl) 

0(Ah x ) 

0(Ah x ) 

0(Ah x ) 

NC-WLSQ 

o(hl) 

o(hl) 

o(h x ) 

0(h x ) 

0(Ah x ) 

CC-SA 

0(hl) 

o(hl) 

0(Ah x ) 

0(Ah x ) 

0(Ah x ) 

CC-NN 

0(hl) 

o(hl) 

0(Ah x ) 

0(Ah x ) 

0(Ah x ) 

CC-FA-unweighted 

0(hl) 

o(hl) 

0(Ah x ) 

0(Ah x ) 

0(Ah x ) 

CC-FA-weighted 

o(hl) 

o(hl) 

o(h x ) 

0(Ah x ) 

0(Ah x ) 

CC-NA 

0(hl) 

o(h x ) 

0(Ah x ) 

0(Ah x ) 

0(Ah x ) 


C. Convergence of discretization errors 




(a) Grids of type ( III P ). 


(b) Grids of type (IV P ) 


Figure 12. Convergence of discretization errors for solution U = sin (nx + 2 Try) on stretched grids of types ( 1 1 1 p j and (IV v ). 

A poor gradient reconstruction accuracy, however, does not necessarily imply large discretization error. Second- 
order accurate solutions have been previously reported l ' 25 on grids with large gradient reconstruction errors. Here, we 
observe similar results for cell-centered and node-centered FVD schemes for constant-coefficient convection. Con- 
vergence histories of the L\ norms of discretization errors for the manufactured solution U = sin (nx + 2ny) on a 
sequence of consistently refined stretched grids of types ( III p ) and ( I V p ) are shown in Figure 12. On grids of type 
(I I I p ), all discretization errors converge with second order. Note that, from the convergence results reported in Sec- 
tion VII (subsection E), discretization-error convergence order for the CC-NA-CLIP scheme is expected to deteriorate 
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to first order on finer grids. Discretization errors of the NC-WLSQ scheme are not shown in Figure 12 because the 
NC-WLSQ scheme does not converge in DCI on grids of types ( III p ) and ( / V p ) . The NC scheme converges with 
first order, as expected. Discretization errors of all cell-centered schemes converge with second order, close to each 
other and to the ideal discretization errors (CC-EG). 


D. Convergence of defect-correction iterations 

The DCI method applied to NC—WLSQ and CC—FA—WLSQ schemes diverges on perturbed stretched grids with 
triangular elements (types (J/ p ), {IIIp), and {IV P ))\ the method converges fast for all schemes on unperturbed grids 
of types (!) — {IV). Somewhat surprisingly, in rectangular geometry, no convergence problems have been detected 
for the CC-NA scheme. Convergence rates of DCI for stable schemes are similar to those observed on isotropic grids 
(Figure 8). Figure 13 shows convergence histories on a 33 x 129 grid of type {IV P ). The asymptotic rates for all 
converging schemes are around 0.5 per iteration. 




(a) Grids of type (1 1 1 p ). 


(b) Grids of type ( IV P ) 


Figure 13. Convergence of L \ -norms of residuals in DCI for FVD schemes with first-order drivers on stretched grids of types ( III P ) and 
(IV P ) with maximum aspect ratio A = 10 3 . 


IX. Grids with curvature and high aspect ratio 

In this section, we discuss accuracy of FVD schemes on grids with large deformations induced by a combination 
of curvature and high aspect ratio. The grid nodes are generated from a cylindrical mapping where (r, 0) denote 
polar coordinates with spacings of h r and hg, respectively; the innermost radius is r = R. The grid aspect ratio is 
defined as the ratio of mesh sizes in the circumferential and the radial directions, A = The mesh deformation is 
characterized by the parameter T: 


t, _ R{l-cos{hg)) _ Rh 2 g _ A hg 

h r ~ 2h r 2 ■ ( ’ 

The following assumptions are made about the range of parameters: R ss 1, A 1, and Th r <C 1, which implies 
that both h r and hg are small. For a given value of A, the parameter I may vary: r > 1 corresponds to meshes 
with large curvature-induced deformation; T <C 1 indicates meshes that are locally (almost) Cartesian. In a mesh 
refinement that keeps A fixed, T = 0{Ahg) asymptotes to zero. This property implies that on fine enough grids with 
fixed curvature and aspect ratio, the discretization error convergence is expected to be the same as on similar grids 
generated on rectangular domains with no curvature. 

We focus on convergence of discretization errors on high-T grids with large curvature-induced deformations. 
Considered manufactured solutions predominantly vary in the radial direction of small mesh spacing. 
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Four basic types of 2D grids are studied in the cylindrical geometry. In distinction from the computational grids 
used in the rectangular geometry, random node perturbation is not applied to high-r cylindrical grids because even 
small perturbations in the circumferential direction may lead to non-physical control volumes. 



(a) Grid of type (III). 



Figure 14. Representative 9 x 33 stretched high-r grids. 

Computational grids are stretched grids with radial extent of 1 < r < 1.2 and angular extent of 20° with a 
fixed maximal aspect ratio A « 1, 100. The grids have four times more nodes in the radial direction than in the 
circumferential direction. The maximal value of parameter T changes approximately from 24 to 1.5. The stretching 
ratio is changing as (3 = 1.25, 1.11, 1.06, 1.03, and 1.01. Representative stretched grids of types {III) and {IV) are 
shown in Figure 14. The tests are performed for the manufactured solution U = sin(57rr). 

A. Approximate mapping method 

Computations and analysis reported earlier 23 ' 25 26 conclude that the unweighted-least-square gradient approximation 
is zeroth order accurate on deformed grids with high F. To improve the accuracy of gradient reconstruction, a least- 
square minimization in a mapped domain is proposed. A general approximate mapping (AM) method based on the 
distance function has been introduced in Ref. [1], 

The AM method applies the LSQ minimization in a local coordinate system, (£, 77), where ?/ is the coordinate 
normal to the boundary and £ is the coordinate tangent to the boundary. The unit vector normal to the boundary, ho, is 
constructed using the distance function, readily available in practical codes, as 

no = ( r o — r o)/ |r 0 — TqI , (15) 

where the position of the control-volume center is denoted ro and the position of the closest point on the boundary is 
denoted rjj. The unit vector tangent to the boundary is denoted as to- 

For constructing the least-square minimization at a control-volume with the center ro, the local coordinates of a 
stencil point r, are defined as 


ii = (r i — r 0 ) -t 0 , 


(16) 


Vi = {Si - So)- (17) 

where Si denotes the distance function of location r,. Thus the ^-coordinate corresponds to the distance from the 
boundary and the ^-coordinate is the projection onto the surface. The least-square minimization yields gradients in the 
(£, 77 ) directions or, equivalently, through a coordinate rotation, in the ( x , y) Cartesian directions. 


18 of 24 


American Institute of Aeronautics and Astronautics 




The left and right states at a control-volume face location, say r f, are reconstructed using gradients in the (£, rj) 
directions along with constructed coordinates 


£/ = ( r / - r o) • t 0 , 

(18) 

Vf = (sf - s 0 ). 

(19) 


The coordinate s / should be an accurate approximation to the distance function from the actual surface, reconstructed 
from points on the actual surface and not from the distance function computed at the interface location. A posible 
approximation is 


s f = (s° f +s})/2, (20) 

where, for node-centered schemes, s° and .s j- correspond to the distance function of the two nodes defining the edge, 
and, for cell-centered schemes, and .s j- correspond to the distance function of the two cell centers adjacent to the 
face. For cell-centered schemes, direct reconstruction using Cartesian coordinate gradients is also possible, yielding 
identical results for grids constructed using advancing-layer techniques. As yet, the AM method has been applied only 
to the cell-centered schemes. 

B. Accuracy of gradient approximation 

The accuracy of gradients reconstructed in the global Cartesian coordinate system for the manufactured solution U = 
sin(57rr) on high-T grids of types (I) — (IV) is summarized in Table 6. Convergence of the maximum gradient errors 
over all control volumes is tabulated. 

Only schemes using the WLSQ method are capable of accurate gradient reconstruction on irregular high-T grids. 
The NC-WLSQ scheme reconstructs accurate gradients on deformed grids of all types. All other schemes show large 
0(1) errors on mixed-element grids of type (IV) with T 1. On grids of type (III), the CC-FA-WLSQ also provides 
accuracy for gradient reconstruction. Schemes using unweighted least-square gradient reconstruction produce large 
gradient errors even on regular grids. 

Table 6. High-F grids: relative errors of gradient reconstruction in global Cartesian coordinates. 



(1) 

(II) 

(III) 

(IV) 

NC 

0(1) 

0(1) 

0(1) 

0(1) 

NC-WLSQ 

0(hj) 

0(hj) 

0(hg) 

0(hg) 

CC-SA 

0(1) 

0(1) 

0(1) 

0(1) 

CC-NN 

0(1) 

0(1) 

0(1) 

0(1) 

CC-FA 

0(1) 

0(1) 

0(1) 

0(1) 

CC-FA-WLSQ 

0(hl) 

0(h e ) 

O(hff) 

0(1) 

CC-NA-CL1P 

0(hg) 

0(hg) 

0(1) 

0(1) 

CC-NA 

0(h%) 

0(hg) 

0(1) 

0(1) 


Gradient accuracy is dramatically improved with the AM method. Table 7 shows accuracy orders for gradients 
reconstructed with cell-centered least-square methods in the local coordinates. All tested schemes provide accurate 
gradients on grids of all types. For illustration. Figure 15 shows relative accuracy of gradients reconstructed on grids 
of type (IV). Note that the CC-NA scheme produces very large gradient errors. This behavior can be explained by 
possible node averaging degeneration on high-T mixed-element grids. On these grids, there are topologies where the 
node solution is averaged from four neighboring cells. The four cell centers involved in such averaging may be located 
on a straight line, thus leading to degeneration. 
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Table 7. High-F grids: relative errors of gradient reconstruction in local AM coordinates. 



(I) 

(ii) 

(III) 

(IV) 

CC-SA 

0{h%) 

O(hg) 

0(h e ) 

0(hg) 

CC-NN 

o(hl) 

0(hg) 

0(he) 

0(hg) 

CC-FA 

o(h%) 

O(h 0 ) 

0(hg) 

0(hg) 

CC-FA- WLSQ 

o(h 2 e ) 

O(h 0 ) 

0(hg) 

0(hg) 



(a) Cartesian coordinates 


(b) Approximate mapping 


Figure 15. Convergence of relative gradient errors for FVD schemes on high-r stretched grids of type ( TV) with maximum aspect ratio 

A = 1.100. 


C. Discretization error convergence 

Convergence of Li -norms of discretization errors of FVD schemes with and without approximate mapping is shown 
in Figure 16. Discretization errors of the NC-WLSQ scheme in Figure 16(a) are shown only for grids with relatively 
low F; on grids with higher F, DCI do not converge. With the exception of the CC-NA scheme on high-F grids of 
type (IV), all other schemes show second-order convergence and very similar discretization errors. Large erratic dis- 
cretization errors of the CC-NA scheme are probably caused by degeneration of the node-averaging stencil mentioned 
in the previous section. This explanation is supported by the evidence of accurate solutions obtained with the CC-NA 
scheme on low-r grids and on triangular grids of type (III), where such degeneration is impossible. On grids of 
the same size, the discretization errors of schemes using the AM method show less variation and are smaller than the 
errors of the corresponding schemes that do not use the AM method. The level of discretization errors obtained by 
the schemes with 0(1) error in the gradient reconstruction is not much different from the discretization error level 
obtained by the schemes with either the AM method (and first-order accurate gradients) or the exact gradient. 

D. Convergence of defect-correction iterations 

Convergence rates of DCI on irregular high-F grids are shown in Figure 17. The DCI method diverges for the CC-NA 
scheme on grids of both types and for the NC-WLSQ scheme on grids of type (III)', on grids of type (IV), the 
NC-WLSQ scheme slowly converges. Note that for all schemes, beside the CC-SA and CC-FA schemes, convergence 
rates of DCI are slower than the rates on perturbed non-curved grids of similar sizes (compare Figures 13 and 17). 
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(a) Grids of type ( III ); Cartesian coordinates 


(b) Grids of type ( IV ) ; Cartesian coordinates 



(c) Grids of type (IV); approximate mapping 


Figure 16. Convergence of L\ -norms of discretization errors of FVD schemes on high-r stretched grids with maximum aspect ratio 

A = 1,100. 


X. Conclusions 

Two node-centered and six cell-centered schemes have been compared for finite-volume discretization of a constant- 
coefficient convection equation as a model of the inviscid flow terms. The cell-centered nearest-neighbor (CC-NN) 
scheme has the lowest complexity; in particular, its stencil involves the least number of neighbors. A version of the 
scheme that involves smart augmentation of the least-square stencil (CC-SA) has only marginal complexity increase. 
All other schemes have larger complexity; the complexity of node-centered (NC) schemes are somewhat lower than 
complexity of cell-centered node-averaging (CC-NA) and full-augmentation (CC-FA) schemes. Defect-correction 
iterations (DCI) with a first-order driver is typically used for solutions of second-order finite-volume discretization 
(FVD) schemes. Convergence of DCI is an important consideration. The CC-NN and CC-SA schemes are promising 
as candidates to be iterated with full second-order linearization. 

Comparisons of accuracy and convergence rates of DCI have been made for two classes of tests: the first class is 
representative of adaptive-grid simulations and involves irregular grids with discontinuous metrics; the second class is 
representative of high-Reynolds number turbulent flow simulations over a curved body. All tests have been performed 
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(a) Grid of type (III). 


(b) Grid of type (IV) 


Figure 17. Convergence of Li -norms of residuals in DCI for FVD schemes with first-order drivers on high-r stretched grids with maximum 
aspect ratio A = 1 , 100. 


for smooth manufactured solutions. 

For the tests of the first class performed in rectangular geometries on consistently refined grids with discontinuous 
metrics, the following observations have been made: 

(1) Discretization errors of second-order schemes are quantitatively similar on grids with the same number of de- 
grees of freedom. The demonstrated convergence of discretization errors closely approaches an “ideal” second- 
order convergence on given grids exhibited by the cell-centered scheme with exact gradients. 

(2) As expected, the NC discretization errors converge with second order on triangular and regular quadrilateral 
grids and with first order on mixed-element (types (IV) and ( I V p ) ) and perturbed quadrilateral (type (Ip)) 
grids. 

(3) Discretization errors of five of the six cell-centered schemes, CC-NN, CC-SA, CC-FA, CC-FA-WLSQ, and 
CC-NA, converge with second order on all tested grids. 

(4) The CC-NA scheme with clipping (CC-NA-CLIP) fails to approximate gradients and degrades solution accuracy 
to first order. The deterioration of solution accuracy is observed on very fine grids with an increased percentage 
of clipped nodes. On coarser grids, the accuracy of the clipped solutions is similar to the accuracy of other 
second-order schemes. 

(5) All schemes may produce 0(Ah x ) large relative errors in gradient reconstruction on perturbed grids of types 
(I p ) — (IV P )\ here A is the grid aspect ratio and h x is the larger mesh spacing. 

(6) As expected, truncation error convergence order is typically one order lower than the convergence order of 
corresponding discretization errors. 

(7) The DCI method for FVD schemes employing weighted least-square gradient reconstruction (CC-FA-WLSQ 
and NC-WLSQ) diverges on perturbed stretched grids. DCI convergence rates for all other schemes, including 
CC-NN and CC-SA, are very fast, while slightly grid dependent; the asymptotic convergence rate is typically 
better than 0.5 per iteration. 

(8) As a recommendation for computations in geometries with no curvature, cell centered CC-NN and CC-SA 
schemes offer the best options of the lowest complexity and second-order discretization errors. 

The tests of the second class have been performed on consistently refined stretched grids generated around a curved 
body, typical of those generated by the method of advancing layers. The range of grid parameters has been chosen 
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to enforce significant curvature-induced grid deformations, characterized by the parameter 1’. All tests have been 
performed for a manufactured solution smoothly varying in the radial direction. 

(1) The discretization errors converge with second order and are small (approaching “ideal” second-order errors) 
for the CC-NN, CC-S A, and CC-FA schemes on all grids and for NC schemes on triangular grids. The errors are 
similar on grids with the same number of degrees of freedom. The discretization errors of the CC-NA scheme 
without clipping do not converge on irregular high-T grids. 

(2) The CC-NN, CC-SA, and CC-FA schemes with least-square gradient reconstruction performed in local ap- 
proximate mapping coordinates provide accurate gradients on all grids. Approximate mapping accounts for 
the global curvature and relies on the distance function that is typically available in practical computations. 
With least-square gradient reconstruction performed in global Cartesian coordinates that do not account for 
global curvature, only the NC-WLSQ scheme provides accurate gradients on all grids; all other schemes fail for 
mixed-element grids of type (IV), generating 0(1) errors in gradient reconstruction. On grids of type (III), the 
only cell-centered scheme with accurate gradient is CC-FA- WLSQ scheme. Note that unweighted least-square 
schemes fail to approximate gradients even on regular grids of types (I) and (II). CC-NA schemes provide 
accurate gradients on regular grids, but exhibit poor gradient accuracy on irregular grids, even with approximate 
mapping. 

(3) The DCI method for the CC-NA scheme without clipping diverges; for the NC-WLSQ scheme, the method 
either diverges or converges very slowly. Convergence rates of DCI for the CC-SA and CC-FA schemes are fast 
and almost grid independent; the average convergence rate is better than 0.5 per iteration. The DCI convergence 
rates for other schemes are slower. 

(4) As a recommendation for computations in curved geometries, the best option is the CC-SA scheme that offers 
low complexity, second-order discretization errors, and fast convergence of DCI. The CC-NN is a promising 
candidate to be iterated with full second-order linearization. The approximate mapping provides uniform accu- 
racy for gradient reconstruction. 
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